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1.  INTRODUCTION 


The  classic  problem  in  signal  processing  is  to  identify  a  low  strength  signal  in  a  noisy 
background.  Even  more  difficult  problems  are  to  identify  several  signals  in  the  noisy  back¬ 
ground,  to  determine  the  direction  of  the  sources,  and  to  track  one  or  more  of  these  multiple 
sources.  In  the  ne.xt  section  we  will  give  an  overview  of  the  techniques  that  we  have  devel¬ 
oped  for  nonlinear  signal  identification.  In  the  third  section  we  will  provide  the  details.  In 
the  fourth  section  we  will  give  an  overview  of  our  nonlinear  tracking  methods,  and  in  the 
fifth  section  we  will  give  the  details  of  the  design  of  an  array  suitable  for  use  with  broadband 
signals.  In  the  sixth  section  we  will  present  the  details  of  our  nonlinear  tracking  algorithm 
for  following  a  maneuvering  vehicle. 

Classical,  linear  theory  has  becii  u^ed  extensively  in  solving  these  problems.  However, 
if  the  signals  have  a  broad  frequency  spectrum  with  no  readily  identifiable  characteristics, 
standard  spectral  methods  are  inadequate.  In  Phase  I,  we  have  built  upon  some  results  of 
nonlinear  dynamics  to  develop  nonlinear  algorithms  to  solve  these  problems,  namely; 

•  classify  sources  of  multiple  signals  in  a  noisy  environment; 

•  determine  the  source  direction  if  the  receiver  is  an  array; 

•  perform  tracking  of  a  maneuvering  vehicle. 

Imagine  there  are  several  sources  of  signals  present  in  a  noisy  environment,  and  it  is 
desired  to  classify  the  sources.  Or  suppose  that  you  have  been  busy  listening  in  a  noisy 
environment  and  a  new  signal  source  suddenly  appears  which  you  wish  to  quickly  identify. 
These  are  classic  problems  for  which  we  developed  a  new  generalized  method  in  our  Phase-I 
effort.  Our  method,  which  is  based  upon  our  research  in  nonlinear  dynamics,  specifically 
addresses  signals  that  have  a  broad  frequency  spectrum  with  no  easily  identifiable  features 
so  that  classical  identification  methods  are  ineffective.  Our  method  is,  to  the  best  of  our 
knowledge,  completely  new,  and  is  based  upon  standard  results  in  nonlinear  dynamics  and 
classical  probability  theory.  It  is  particularly  well  suited  to  signals  associated  with  low¬ 
dimensional  chaotic  dynamical  systems,  but  it  is  by  no  means  restricted  to  them. 

VVe  now  give  a  brief  overview  of  the  procedures  and  assumptions  involved  in  performing 
the  classification.  We  assume  that  samples  of  signals  have  been  previously  acquired  and  are 
available  in  a  library.  The  signal  representing  the  noise  background  may  be  obtained  while 
one  is  listening  for  other,  more  interesting  signals,  and  it  is  treated  the  same  as  any  other 
signal.  Certain  geometric  structures  (probability  densities  and  their  characteristic  functions) 
associated  with  each  signal  then  are  constructed  and  a  comparison  is  made  with  the  same 
structure  created  from  the  newly  received  signal.  The  way  in  which  this  comparison  is  made 
can  have  a  significant  effect  on  performance.  So  far  we  have  relied  on  either  a  standard 
statistic  for  comparing  two  probability  densities,  or  upon  a  simpler  square  of  the  difference 
of  the  densities  summed  over  all  signal  strengths. 

The  geometric  structures  are  defined  so  that  the  structure  from  the  sum  of  two  signals 
factorizes  into  the  product  of  structures,  one  from  each  signal.  VVe  then  compare  with  the 
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corresponding  product  of  structures  from  our  library  and  decide  whether  we  have  a  correct 
identification.  This  factorization  property  is  crucial  as  it  permits  a  large  library  of  signals  to 
be  compared  with  the  incoming  signal  in  essentially  real  time.  The  result  is  a  determination 
that  the  received  signal  has  a  component  that  was  made  by  a  source  of  a  unique  type 
such  as,  for  example,  a  particular  helicopter  model  or  a  particular  submarine  model.  For 
spectrally  wide-band  chaotic  signals  it  is  usually  not  possible  to  make  such  an  identification 
by  comparing  either  signals  or  their  spectra. 

One  can  imagine  using  this  kind  of  processing  in  different  ways.  One  scenario  involves 
listening  at  a  particular  site  with  a  nearly  stationary  (in  time)  noise  background.  Suddenly 
a  new  noise  is  present,  and  it  is  desired  to  identify  the  source.  In  this  case  we  need  only 
identify  the  new  signal  as  we  can  include  the  observed  background  noise  in  our  library.  The 
searching  of  the  library  is  then  very  rapid  and  an  identification  may  be  quickly  performed. 
In  other  situations  one  may  be  confronted  with  a  new  environment  and  wish  to  identify  all 
ot  the  sources  present. 

In  some  situations  the  receiver  is  an  array,  and  the  signals  from  the  individual  elements 
can  be  separately  phased  (independent  time  lags)  before  combining  the  signals.  In  this  situ¬ 
ation  it  is  possible  to  identify  the  direction  of  a  particular  source  by  varying  the  phase  of  the 
array.  Because  the  signal  is  broadband,  we  are  not  looking  for  a  null  to  determine  the  direc¬ 
tion,  but  rather  changes  in  the  constructed  geometric  structures  as  the  phases  of  the  array 
elements  are  changed.  To  the  best  of  our  knowledge,  this  is  a  new  concept.  The  densities  con¬ 
structed  from  noise-like  signals  were  observed  to  have  relatively  little  structure  as  the  phase 
is  changed,  but  the  densities  from  other  signals  showed  significant  variation  as  a  function  of 
the  array  phase.  It  is  also  possible  to  monitor  changes  in  the  distance  of  the  source  as  the 
amplitude  of  the  received  signal  varies.  new  algorithm  that  we  have  developed  determines 
the  velocity  of  the  source  by  a  nonlinear  “doppler  shift”  calculation.  This  determination  is 
sensitive  to  the  background  noise;  consequently,  its  practical  implementation  will  have  to  be 
delayed  until  the  nonlinear  signal-processing  algorithms  are  refined. 

In  conjunction  with  the  development  of  our  nonlinear  identification  algorithms,  we  have 
also  derived  a  nonlinear  tracking  algorithm.  It  also  turns  out  to  be  suitable  for  tracking 
maneuvering  vehicles  in  a  conventional  signal-processing  environment. 

In  Phase  I  we  performed  the  following  tasks: 

•  VVe  constructed  twenty  signals,  including  several  types  of  broadband  noise,  outputs  of 
chaotic  electronic  circuits,  and  samples  of  time  series  generated  by  systems  of  ordinary 
differential  equations  with  chaotic  solutions.  Samples  of  the  signals  were  placed  in  a 
“library.”  Additional  signals  from  the  same  systems  were  generated  and  constitute 
the  “test”  set.  The  signals  in  the  test  set  were  not  identical  to  any  in  the  library, 
but  they  were  from  the  same  dynamical  systems.  In  addition,  we  forced  some  of 
the  noise  signals  to  have  the  same  spectrum  as  some  of  the  chaotic  signals.  This 
was  accomplished  by  Fourier  transforming  the  signal  and  randomizing  the  phases  of 
the  Fourier  components  and  then  inverting  the  Fourier  transform.  This  insures  that 
standard  spectral-identification  methods  could  not  separate  the  sources. 


Two  signals  from  tlie  test  set  were  added  together  with  randomly  selected  amplitudes 
to  form  a  simulated  received  signal.  We  then  attempted  to  identify  both  the  amplitudes 
and  the  type  of  signals  by  comparing  certain  geometric  structures  constructed  from  the 
signal  and  from  the  library.  Provided  the  amplitudes  of  the  test  signals  were  not  too 
different,  we  were  able  to  identify  both  signals  in  all  cases.  The  range  of  amplitudes 
permitted  depended  upon  the  similarity  of  the  signals,  with  a  factor  of  roughly  ten 
being  a  limit  for  the  modest  amount  of  processing  that  we  did.  We  believe  that  as  the 
algorithms  are  refined,  substantial  improvements  in  the  range  of  amplitudes  for  which 
signals  can  be  separated  will  be  possible. 

•  We  also  combined  three  test  signals  together  for  a  few  e.xamples  and  were  able  to 
classify  the  components. 

•  We  next  considered  an  example  where  the  receiver  was  a  multi-component  array.  The 
incoming  signals  from  the  receivers  in  the  array  were  added  with  a  phase  difference 
which  was  chosen  to  select  a  particular  direction.  Because  the  signals  are  broadband, 
there  is  no  "nuU'’  as  a  function  of  the  relative  phases  of  the  array  elements.  However, 
the  constructed  densities  have  a  strong  dependence  on  the  phase,  which  can  therefore 
serve  both  as  a  direction  identifier  and  as  an  additional  classification  mechanism.  In 
examples  of  arrays  that  were  typical  of  ones  that  could  be  deployed  in  the  ocean, 
we  were  able  to  determine  the  direction  and  identify  the  source  of  signals  that  were 
20db  weaker  than  other  signals  with  similar  spectral  characteristics  with  only  a  modest 
amount  of  processing. 

The  techniques  that  we  developed  are  new,  and  there  are  a  wide  variety  of  possible 
refinements.  Because  the  processing  is  nonlinear,  there  are  more  variations  available  than 
for  conventional  linear  signal  processing,  even  when  we  restrict  the  investigation  to  the  types 
of  algorithms  that  we  propose.  There  are  no  real  differences  in  the  principles  involved  in 
the  various  implementations,  but  there  are  many  different  ways  to  optimize  the  probability 
of  correct  classification.  Unfortunately  it  is  very  difficult  to  prove  general  theorems  about 
what  is  optimum  for  nonlinear  systems.  Thus  we  expect  that,  for  the  next  few  years  at  least, 
the  non-linear  signal  processing  community  will  have  to  rely  to  a  large  extent  on  empirical 
evidence. 

2.  OVERVIEW  OF  OUR  SIGNAL  CLASSIFICATION  METHOD 

We  now  discuss  some  of  the  important  variables  and  some  of  the  freedom  available  in 
the  details  of  our  algorithms.  Usually  our  signal  is  a  time  series  from  a  single  detector. 
In  nonlinear  dynamics  we  are  often  interested  in  quantities  that  might  have  been  observed 
if  more  detectors  were  used.  In  order  to  get  some  additional  information,  we  construct 
additional  surrogate  signals.  There  is  a  certain  arbitrariness  in  this  procedure  which  leads  to 
a  wealth  of  simple  variations  in  processing  algorithms.  For  example,  the  most  common  way 
of  constructing  a  second  signal  is  to  start  with  the  original  signal,  S{t),  and  construct  another 
signal,  Si{t)  =  S{t  +  T)  where  T  is  a  time  delay.  The  two  signals,  S{t)  and  Si{i)  can  be 
used  to  construct  a  two-dimensional  phase  space.  Other  signals  can  be  constructed  by  using 
other  time  lags,  such  as  S-iU)  =  Sil  +  2T).  For  dynamical  systems  with  low  dimensional 


attractoi-s,  it  is  possible  to  construct  a  phase  space  that  luis  geometric  properties  related  to 
the  actual  phase-space  geometry.  This  tcchniciueof  phase-space  reconstruction  has  been  used 
extensively  to  obtain  information  about  dynamical  systems.  It  underlies  our  procedures,  but 
we  have  found  that  it  can  be  extended  in  several  ways  for  purposes  of  signal  recognition. 

In  most  nonlinear  dynamics  applications,  one  considers  a  fixed  time  delay,  T,  and  for 
those  applications  nothing  is  gained  by  varying  the  time  delay,  except  perhaps  one  choice 
works  slightly  better  than  another.  For  our  applications  it  is  useful  to  consider  several 
time  delays  simultaneously.  Of  course,  this  is  not  a  new  concept  in  linear  signal  processing 
since  the  auto-correlation  function  considers  all  time  delays.  However,  this  idea  has  not 
previously  been  used  for  examining  phase-space  density  functions.  We  found  that  making 
separate  comparisons  for  several  different  time  lags  improved  the  decision-making  algorithm. 
•As  one  example,  an  alternative  way  of  constructing  a  second  signal  is  to  consider  the  signal 
Ssum{t)  =  S{t)  +  S{t  +  T).  Although  this  is  an  apparently  trivial  change,  it  turns  out  to  be  a 
crucial  ingredient  when  the  receiver  is  an  array,  and  it  provides  a  significant  improvement  in 
our  classification  algorithm  for  single  receivers.  The  density  constructed  from  such  a  signal 
combination  is  a  slice  through  a  two-dimensional  density  constructed  from  .?i  and  So-  but 
it  usually  has  more  structure  than  the  density  constructed  from  .Fj  alone  and  is  hence  more 
useful. 

The  point  of  the  above  two  paragraphs  is  that  there  are  a  variety  of  simple  modifications 
that  can  be  made  to  our  basic  method  that  will  improve  the  results.  .As  we  have  indicated 
above,  there  are  several  variations  in  the  numbers  and  types  of  phase-space  densities  that 
can  be  used  for  identification.  The  principal  ones  that  we  have  tested  are  as  follows: 

•  Two-dimensional  phase  space  with  a  single,  fixed  delay  time: 

•  Two-dimensional  phase  space  with  one  or  two  additional  delay  times; 

•  One-dimensional  phase  space  using  Ssum{i)  many  delay  times; 

•  Three-dimensional  sampled  phase  space  with  a  fixed  delay  time. 


Most  of  our  tests  have  been  performed  using  a  two-dimensional  phase  space  with  a  fixed 
delay  time,  T,  and  it  is  within  that  framework  that  we  declare  Phase  I  to  be  a  success. 
In  order  to  improve  the  discrimination,  we  first  considered  a  second  time  lag  and  that  was 
sufficient  to  remove  some  accidental  ambiguities.  It  is  po.ssible,  and  in  principle  desirable,  to 
work  in  a  three-dimensional  or  higher  phase  space.  Unfortunately,  the  computational  time 
increases  dramatically.  Rather  than  reconstruct  the  entire  three-dimensional  density,  it  is 
possible  to  sample  randomly  the  density  (actually  the  Fourier  transform  of  the  density)  at 
many  points.  This  allows  one  to  enjoy  many  of  the  benefits  of  a  full  reconstruction  without 
the  full  computational  effort.  We  performed  some  limited  testing  using  this  procedure  in 
Phase  I  and  it  appc.'ars  to  be  successful.  If  several  signals  are  present,  this  technique  might 
provide  a  better  signal  to  noise  ratio  for  the  same  computational  effort.  We  have  not  yet  made 
comparisons  with  the  use  of  two-dimensional  densities  in  terms  of  efficiency  and  accuracy. 
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We  also  did  some  calculations  with  one-dimensional  densities  using  the  signal  = 

S{t)  -r  S{t  -h  T)  described  above.  The  disadvantages  ot  a  density  constructed  from  just  S{t) 
is  that  it  is  dilhcult  to  obtain  good  discrimination  because  of  the  relativ'ely  few  features 
available  for  discrimination,  and  that  drawback  is  related  to  the  relatively  little  information 
contained  in  that  density.  However,  adding  a  time-lagged  signal  increases  the  structure 
and  considering  many  different  time  lags  greatly  increases  the  identification  capability.  A 
significant  advantage  of  one-dimensional  densities  is  that  relatively  little  computer  resources 
are  required.  We  have  tested  this  concept  on  a  few  signals.  .At  the  very  least  it  is  successful 
enough  to  serve  as  a  quick  “rejection  filter,”  thus  reducing  the  number  of  cases  requiring 
more  extensive  examination. 

.As  we  mentioned  above,  these  combinations  and  their  generalizations  also  pla}'  an  im¬ 
portant  role  when  the  receiver  is  an  array.  In  the  plane-waves  approximation,  the  best  one 
can  hope  to  do  is  to  translate  the  phase  information  into  a  bearing  angle.  The  problem 
is  thus  to  design  an  array  whose  beam  pattern  has  little  or  no  frequency  dependence.  In 
particular,  we  must  make  sure  that  the  main  lobe  width,  peak  response,  sidelobe  level,  and 
distance  between  the  main  lobe  and  the  sidelobe  plateau  (which  takes  the  place,  in  our 
approach,  of  the  grating  sidelobes  which  appear  genericalK'  in  the  case  of  single  frequency 
arrays)  all  remain  constant  across  the  design  frequency  band.  To  do  this,  we  use  the  Poisson 
summation  formula  and  the  method  of  stationary  phase  to  relate  beam-pattern  properties 
to  array  characteristics.  Using  this  relation,  we  then  translate  all  the  listed  beam-pattern 
requirements  into  practical  requirements  on  sensor  location  and  array  amplitude  shading. 
We  end  up  with  a  beamforming  approach  that  is  significantly  more  efficient  than  a  classical 
uniformly  spaced  array.  Specifically,  while  uniform  spacing  produces  an  array  in  which  the 
total  number  of  elements  would  be  proportional  to  the  ratio  ot  the  highest  to  lowest  design 
frequencies,  the  method  we  have  u.sed  allows  one  to  achieve  all  the  beam-pattern  objectives 
with  a  total  number  of  elements  that  grows  like  the  logarithm  oi  the  ratio  of  highest  to  lowest 
frequencies. 

In  summary,  in  Phase  I  we  ha\e  developed  new,  nonlinear  algorithms  for  detecting  one 
or  more  spectrally  broadband  signals  in  a  noisy  background.  The  primary  goal  of  verifying 
that  signals  can  be  identified/classified  by  the  proposed  technique  was  achieved.  In  addition, 
the  fact  that  there  are  several  variations  of  constructing  densities  that  have  been  successful 
in  identification  seem  to  enhance  greatly  the  ideas’  usefulness.  The  algorithms  use  a  new 
and  innovative  approach  that  can  be  used  in  near  real  time.  We  have  successfully  used  the 
algorithms  to  classify  sources  and.  with  an  array  of  receivers,  we  can  identify  the  direction 
as  well  as  decrease  the  contribution  of  other  sources. 

3.  SIGNAL  CLASSIFICATION 

The  crucial  ingredients  for  a  correct  classification  of  signals  are  their  probability  densi¬ 
ties.  First,  we  explain  two  ways  of  constructing  the  densities,  and  then  we  will  show  how 
to  identify  the  components  of  a  signal  made  up  of  several  pieces.  Examples  will  be  given 
showing  how  the  method  works,  and  some  important  refinements  will  be  presented. 

The  first  step  is  to  normalize  all  signals.  The  mean  is  removed  and  the  root-mean- 
sejuared  signal  is  normalized  to  unity.  Il  Is  essential  to  have  a  uniform  normalization  since 
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we  have  no  a  priori  knowledge  of  tlic  relative  amplitudes  of  the  components  of  a  signal.  The 
normalization  and  mean  may  be  useful  themselves,  such  as  for  determining  the  change  in 
distance  of  a  moving  vehicle,  but  we  don’t  use  them  directly  for  most  signal  identification. 

The  simplest  density  is  ju.st  a  histogram  of  the  amplitude  of  the  incoming  signal.  If  S{t) 
is  sampled  and  has  the  values  .S'l,  S-y.  S:i . . .  for  a  total  of  N  values,  we  just  count  the  numljer 
of  times  .9  is  between  x  and  ,r  +  r/.r  and  call  that  N  ■  p[x)dx.  Formally, 

T 

p{x)  =  \\m  6{x  -  S{t))dt  (1) 

wliere  d[x)  is  the  delta  lunction.  This  density  usually  doesn't  have  enough  structure  to  be 
a  uselul  identifier,  although  with  a  modification  to  be  described  below  it  can  be  a  useful 
tool.  .A  second  density  now  two  dimensional,  can  be  constructed  by  introducing  a  second 
coordinate  S{t  4-  7’).  Then 
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The  practical  ap[)!ication  of  these  formulae  requires  replacing  the  delta  function  by  a 
function  of  finite  width.  VVe  have  chosen  to  u.se 
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Other  substitutions  are  pos.sible. 

The  interpretation  of  p{xxij)  was  discussed  in  the  previous  section.  .A  third  coordinate 
can  be  introduced  and  a  three-dimensional  density  defined.  The  procedure  could,  in  principle, 
be  continued  to  whateve’-  dimension  one  wished.  There  are  several  reasons  for  truncating  at 
a  low  dimension: 

•  Computational  requirements  become  prohibitive  in  high  dimensions  (greater  than  three) 

•  The  length  of  the  signal  may  not  be  long  enough  to  justify  many  dimensions,  i.e.. 
regions  of  phase  space  will  be  empty  or  inadequately  sampled. 

•  For  low-dimensional  dynamical  systems  the  higher  dimensional  probability  densities 
are  certainties  in  terms  of  lower  dimensional  ones. 

•  .A  density  in  a  low  dimension  is  a  projection  of  a  density  in  a  higher  dimension.  The  low¬ 
dimensional  density  may  contain  sufficient  features  to  serve  as  a  classifier  for  practical 
purposes. 

It  is  possible  to  sample  a  high  dimension  without  computing  the  entire  density.  This 
procedure  will  be  discussed  below.  For  dynamical  systems  with  low-dimensional  attractors. 
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strange  or  otherwise,  there  are  lln'orems  that  relate  the  luunbi'r  ot  time  lags  needed  to  the 
dimensionality  of  the  system.  Con.secinently,  if  only  low-dimensional  dynamical  systems  were 
the  sources  of  the  signal,  we  need  consider  only  lower  dimensional  probability  distributions. 
Noise  tends  to  be  very  high  dimensional,  but  it  usually  has  little  structure  anvwvay  so  not 
much  is  lost  by  projecting  it  onto  a  lower  dimensional  space. 

We  have  electeil  to  do  our  calculalions  mostly  :n  a  two  dimensional  space.  .A.S  wc  have 
noted,  the  various  time  lags  are  parameters  ol  the  density  lunctions.  To  work  with  only  one 
time  lag,  2’.  is  (.'cjiii\'aleut  to  considering  the  auto-correlation  iunction 

C(T)  =<  Sit)Sit  +  T)  >  (-1) 

at  only  one  value  of  T.  Of  cours('  one  value  ot  C{T)  provides  very  little  useful  information. 
For  e.xample.  the  (.h'lei  ininatiou  of  the  spectrum  requires  a  witle  range  ol  values  of  7  (in  prin¬ 
ciple  all  T  values,  in  practice  a  limited  range).  Thus  even  when  we  restrict  our  calculations 
to  two  dimensions  we  can  construct  a  continuum  of  densities  ilepending  on  7  .  In  jjractice 
we  have  samphxl  a  range  of  T  values. 

The  ne.xt  step  in  the  procedure  is  to  understand  the  relationship  between  a  density  of  a 
signal  and  the  densit\'  cf  its  components.  To  this  end  we  define  the  characteristic  Iunction 
of  the  density  which  is  just  the  Fourier  transform, 

/;(/.•)=  I  c  .r  )(/.r  =  '*<//,  ( •') ) 

This  has  the  property  that  if  S{t)  =  -?i(0  +  S2{t)  and 

Mk)  =  j.  [ Ci) 

then  p{k)  =  pi{k)f)-2{k).  This  is  a  well-known  property  of  characteristic  functions  [5]  and 
just  requires  that  .S’l  and  S-2  are  generated  by  independent,  stationary  processes. 

We  emphasize  that  this  relationship  is  crucial  to  the  possibility  ot  signal  classification  if 
multiple  signals  are  present.  If  thrc'e  signals  are  present  the  result  is 

m  ^  (S) 

This  e.xpression  would  be  correct  for  general  signals  .9]  and  .^i,  but  since  we  normalize 
signals,  the  amplitude  must  ent(’r  tlu'  relationship  bolwec'ii  the  ps  in  Ecp  S.  Suppose 


S{t)  =  aSi{t)  +  l3S2it) 


where  all  signals  have  zero  mean  aiul  are  properly  normalized.  Tiien  a’  +  .i"  =  !.  and  tlu; 
relationship  is  actnally 


p(  k]  =  /)  1  ( o  A- )  /)  j  ( 3  A- )  .  (Id ) 

Since  a  and  -i  aren't  known,  this  re(|nires  that  they  be  I’onnd  as  part  of  the  procedur*'. 
The  relationship  is  the  sanu.'  in  high<‘r  dimensions  e.\e(,“pt  that  A’  is  now  a  vector. 

1  he  density  can  Ije  ciinst  rncti-d  either  by  const  met  ing  f){x  )  or  p(A-)  from  tlu’  above 
delinitions.  It  is  nsnally  more  elliciv-nt  to  calculate  p(A-)  on  a  mi'sh  ami  then  interpolat<' 
to  get  /TiWy),  ratlic-r  than  calcul.ite  directly.  .\ote  that  p{.v)  will  always  vanish  for  .r 

outside  some  region,  since'  signals  can't  Inive  arbitrarily  large  amplitudes. 

Onci'  tlu;  ;i[)pr(jpriat  densities  are  const  ructeel.  the  ne.xt  ste[)  is  to  te.'st  h\'i)ot  lu'si's  aluait 
which  signals  make  up  .'t'’.  The  procedure  that  we  have  used  thus  far  to  ]>rovide  a  figure  of 
merit  is  to  test  the  (piantit}' 


Ifrror  =  ^  |  /7(A-)  -  pi(oA-)/3j(.AA-)  |' 

k 


(II) 


or  the  corresponding  Fourier  transform  as  a  futiction  ol  .r  lor  all  signahs.  .Fj  and  5j.  for  a 
range  of  values  of  a.  If  one  of  the  components  is  known  then  one  need  only  search  for  a 
second  signal.  This  is  the  cas<>  if  a  new  signal  suddenly  appears.  The  background  signal 
up  to  that  tiiiK'  is  5i  and  tlu'  new  signal  is  5-2-  Since  we  have  been  listening  to  5i.  we 
can  construct  tlu’  densitic's  for  it  and  only  candidates  for  5-2  need  to  be  tested.  This  also 
illustrates  that  one  of  the  signals.  .Fj  in  this  case,  can  be  arbitrary  noise.  If  there  are  three 
signals  present,  all  unknown,  one  has  to  test  many  more  hypotheses. 

We  have  found  that  a  s[)<;cial  class  of  one-dimensional  densities,  to  be  described  below, 
can  be  used  as  a  rejection  filter.  This  is  an  important  point  since  they  require  very  little 
com p u t at i on al  resou rces . 

We  now  present  e.\am[)les  to  illustrate  the  method.  In  Figures  1  and  2,  several  e.xamples 
of  signals  are  shown.  Two  of  the  signals  are  from  electronic  circuits  (provided  by  NRL) 
and  the  rest  are  s\  nfhetic  signals  generated  fioni  ordinary  difTerential  equations  with  chaotic 
attractors.  Figures  3  and  4  show  the  power  spectrum  of  the  same  signals. 

The  basic  idea  of  the  scheme  we  propose  for  classifying  signals  is  to  construct  probability 
densities  for  signals  of  interest  and  to  make  comparisons  by  pattern  matching  with  the 
corresponding  density  for  a  signal  that  one  wishes  to  identify.  There  are  an  infinite  number 
of  possible  densities,  but  it  is  sufficient  to  compare  only  a  few.  Before  proceeding  to  the 
technical  details  we  give  a  pictorial  example.  In  Figure  1  there  are  three  density  plots.  The 
first,  (a),  is  tlu'  density  constructed  from  a  signal  to  be  classified  and  the  other  two.  (b) 
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and  (c),  are  [proposed  candidates  taken  from  a  library.  In  tliis  case  it  is  easy  to  see  that  (c) 
matches  (a)  and  that  (b)  do(.'s  not.  .A  numerical  calculation  of  the  differences  confirms  this. 
!  '  :  s  particular  density  has  a  simple  interpretation.  It  is  the  probability  that  a  measurement 
oi  a  signal,  S(t),  at  time,  f,  and  a  sulxsequent  measurement  a  time  T  later  will  have  the  values 
X  and  ij  respectively.  Or  more  properly,  }j)(lxdij  is  the  probability  for  the  signal  S{t)  to 
be  in  the  range  :r  to  x  -f  dx  and  the  signal  S[l  +  T)  to  be  in  the  range  y  to  y  +  dy  .  p{x.y) 
is  tin?  vertical  a.x'is  lafjeled  density.  The  rlensities  will,  of  course,  be  different  for  differ<?nt 
\'alues  ot  T.  1  lu?  sc'cond  moment  of  the  probabilit}'  density.  p{x.y).  is  the  antocorr<?lat ion 
function,  i.i?..  L\T)  =<  xy  >. 

Sncii  a  simple  t:omparison  do<.'s  not  always  i)rovide  a  conclusive  test,  p;;rticn!arly  when 
mult i[:)le  signals  and  noise  are  irrestuit.  In  that  case  one  can  compare  the  pattc'rns  at  Sfweral 
different  values  of  /  or  construct  other  relateil  densit.ies.  There  are  many  finiinitef  densities 
that  can  b('  const  met  ('il  and  it,  mu'^t  be  de!ermin(?d  which  an'  more  useful  for  a  particular 
application. 

■fTiis  e.xamph'  illustrates  sonu'  of  tlu'  Ih'.xibility  available  in  onr  algorithms  for  signal 
discrimination,  .and  it  also  s('r\'es  to  illustrate  sonu' ot  the-  jM'oblems  th;it  need  to  be  soU’cd  in 
orch'i'  to  de\'elop  a  uselid  s\'st('m.  1  lu're  is  llexibility  in  the  choice' of  (.lim('nsion  of  probability 
density,  there'  is  the'  eiucstiejii  ol  .seh.'ctiexi  eaf  the'  lime  lag.  7'.  and  whether  to  use  many  dilfen'iit 
\’alues  of  T  aiul.  il  so.  how  many.  There  is  the  possibility  of  using  an  array  of  receive'rs  to 
enhance  recognition.  .Since?  the'  sources  have' a  broad  freciuency  spectrum,  the?  eiuestion  arisets 
as  to  t  lie  most  e'lfe'cti\’c  way  oi  deploying  an  array  and  what  are  the  trade-oifs.  Otie  also  neeels 
an  automateej  elecision  making  process  for  comparing  e.lensilic.s  aiul  some  way  to  measure' 
confielence.  .All  of  the  sugge.'stions  for  enhancing  performance  were  tested  successfully  in 
Ph.ise  I.  although  no  information  about  optimizing  performance  was  oblaineid. 

In  order  to  test  our  methoel  when  multiple  signals  are  present,  we  form  combinatieins  of 
two  signals  to  give  a  third. 


S{l)  =  aP[t)  +  bQ[t)  (12) 

The  signals  P  and  Q  are  to  be  identified  as  well  as  their  relative  amplitudes.  Sitice  we 
measure  S{t)  we  know  its  level  so  that  we  know  the  value  of  the  sum  of  the  sciuares  of  a  and 
b.  For  the  particular  example  that  we  show  lu're,  a  =  O.TS  and  b  =0.S7. 

We  have  been  using  a  \‘  statistic  as  a  figure  of  merit,  which  is  a  slight  modification  of 
E(p  11  [5]  and  the  appropriate  function  of  the  dc'iisities  of  the  signals  is  comparerl.  In  this 
case  we  const  ruclc'd  a  two-dimensional  (h'lisity  by  using  a  signal  and  a  time-lagged  signal  for 
the  two  com[)onents. 

1  he  criterion  that  we  use  is  exact  in  the  limit  of  an  infinitely  long  signal,  but  in  the 
present  case  a  relatively  short  signal  was  useil.  This  means  that  even  when  the  correct  signals 
with  the  correct  amplitudes  are  comparo'd  there  will  be  some  residual  statistical  error.  In 
the  lower  right  graph  ol  Figure  6,  \'  is  shown  versus  the  relative  strength  parameter,  o. 
Theix  ?  are  two  curves  that  arc  iK'arly  identical  and  a  third  one  that  has  a  larger  \*.  One  of 
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the  two  lower  curves  corresponds  to  the  exact  comparison  and  gives  us  an  estimate  of  the 
statistical  uncertainty  that  places  a  lower  limit  on  x^-  The  second  curve  was  computed  by 
taking  a  different  sample  of  signals  4  and  6  and  comparing  various  combinations  of  them 
with  the  densities  of  the  original  signal.  We  see  that  the  value  of  Q'  is  clearly  about  0.5. 
which  is  correct. 

The  residual  error  is  due  to  the  inherent  lack  of  statistics.  The  upper  curve  corresponds 
to  the  hypothesis  that  one  of  the  signals  was  correct  and  the  other  incorrect.  Note  that 
even  in  this  case  our  algorithm  correctly  identilies  the  amount  of  the  signal  that  was  the 
right  candidate.  We  then  cousitlered  two  different  signals,  neither  of  which  corresponded 
to  the  original  choice,  namely  signals  1  plus  5  (see  lower  left  of  Figure  6).  In  this  case 
was  relatively  Hat  and  well  above  the  levels  shown  on  the  graph,  thus  affirming  that  neither 
signal  was  present  in  the  signal  to  be  identified. 

The  amount  of  data  used  is  of  course  important.  The  computational  effort  is  propor¬ 
tional  to  the  length  of  the  time  the  signal  is  ob.served,  and  there  is  the  logistic  problem  of 
not  being  able  to  collect  an  infinite  amount  of  data  in  a  finite  time.  For  the  examples  shown, 
we  sampled  about  ten  times  per  oscillation  and  included  a  total  of  10,000  points.  For  typical 
machinery  oscillations  of  1,000  rpm,  this  is  about  a  minute's  worth  of  data. 

Depending  upon  the  ty[)es  of  signals  of  interest,  one  may  get  by  with  less  total  samples. 
We  tried  many  different  combinations  of  signals,  ev'en  including  noise  in  the  signal  to  be 
identified,  and  all  of  our  results  were  comparable  to  the  example  presented  above.  We  were 
even  able  to  recognize  that  a  signal  was  made  up  of  two  similar  signals  so  that  it  there  are 
two  sources  of  the  same  signal  we  could  recognize  that  fact. 

It  is  important  to  be  able  to  identify  a  weak  signal  combined  with  a  strong  signal. 
The  limitations  on  this  are  not  known  at  pre.sent  and  will  undoubtedly  be  dependent  upon 
refinements  of  the  algorithm.  Qualitatively  we  have  no  problem  separating  signals  with  a 
ratio  of  energies  of  ten  to  one.  Sometimes  we  can  resolve  power  levels  of  a  hundred  to  one. 
and  we  expect  further  improvement  will  be  po.ssible. 

WT  have  done  some  limited  testing  with  a  signal  made  up  of  three  signals.  In  the  first 
example  one  signal  was  broad-band  noise  (signal  7)  and  the  two  other  were  chaotic  signals. 
The  energy  in  each  signal  was  the  same.  The  correct  combination  was  easily  selectcil  out. 
We  then  tried  a  combination  with  the  energies  in  the  ratio  of  9  to  4  to  1,  with  the  9  being 
the  broad-band  noise.  There  did  not  seem  to  be  any  problem  identifying  the  two  signals 
even  in  the  presence  of  noise.  There  was  one  example  of  an  alternative  signal  giving  as  good 
a  statistic  as  one  of  the  correct  signals.  In  order  to  resolve  this  discrepancy  we  performed 
a  separ.  test  which  completely  resolved  the  uncertainty.  It  is  encouraging  that  almost 
all  incorrect  combinations  were  easily  rejected  in  the  two-dimensional  phase  space,  as  this 
means  that  only  a  few  candidates  might  have  to  be  further  processed,  as  we  did  in  the  above 
example  in  order  to  get  better  rejection  of  incorrect  hypotheses. 

In  order  to  give  the  reader  a  visual  impre.ssion  of  the  densities  we  illustrate  with  a 
combination  of  two  signals,  such  that  the  weak  signal  has  40%  of  the  amplitude  of  the  first 
signal.  We  then  consider  three  candidate  signals  for  the  weak  signal.  Figure  7a  shows  the 
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density  p[x.ij)  ot  the  incoming  coinijosite  signal  with  x  =  S(t)  and  ij  —  5(<  +  4).  T!ie  density 
obtained  Iw  taking  the  convolution  ol  the  two  correct  individual  signal  densities  is  shown  in 
Figure  7b.  Note  that  density  in  Figure  7b  is  smoother  than  Figure  7a.  This  is  due  to  the 
limited  statistics.  It  the  signals  were  infinitely  long,  the  two  densities  woidd  be  ioentical. 
Incorrect  hypotheses  for  the  weak  signal  combined  with  a  correct  guess  for  the  strong  signal 
are  shown  in  Figure's  7c  and  7d.  Note  that  visually  we  can  distinguish  them,  and  a  nun.erical 
test  also  concludes  which  is  the  correct  density. 

\'Ve  now  turn  to  a  promising  method  ol  perlorming  some  classification  via  one-dimensional 
densities.  Define  a  new  signal.  Ii[t.)  =  S{t)  -|-  S{1  T).  Construct  the  density  associated 

with  this  signal. 


1 

P{~r.T)  =  —  8{x- R{t))dt.  (13) 

Jo 

As  usual  we  normalize  R{t)  such  that  nii.s{R^)  =  1.  The  factor  involved  in  the  normal¬ 
ization  is  just  the'  correlation  tunction,  that  is,  rins{R^)  =  2(1  -|-  C'{T))  before  normalizing 
R  but  alter  normalizing  5.  The  density  p{x;T)  has  considerable  structure  and  the  structure 
changes  as  T  changes.  In  Figure  S.  we  present  e.xamples  for  one  particular  signal.  The 
density  is  shown  on  \’crtical  a.xis  and  the  signal  strength,  x,  is  on  the  horizontal  a.xis.  The 
value  of  T  is  shown  at  the  top  ot  each  graph.  For  reference,  a  typical  oscillation  time  scale 
is  around  ten  units.  The  dejisities  for  a  second  signal  are  shown  in  Figure  9.  Note  that  they 
are  quite  distinct,  particularly  it  several  dilferent  T  values  are  examined.  The  density  with 
T  =  0  corresponds  to  the  usual  density.  Note  that  it  is  relatively  featureless.  One  way  of 
constructing  noise  with  a  Ircquency  spectrum  comparable  to  a  signal  is  to  Fourier  transform 
the  signal  and  make  the  phases  random  with  each  one  independent  of  the  others.  Then 
Fourier  transform  back  and  the  resulting  signal  will  have  the  same  spectrum.  The  result 
of  such  an  operation  is  shown  in  Figure  10  a-f.  Notice  that  there  is  no  structure  and  no 
dependence  on  T .  Figure  lOg  is  tor  S{t.)  =  cos(f).  There  is  no  dependence  of  p{x)  on  T  for 
a  simple  sinusoidal  signal. 

We  have  done  some  limited  testing  of  comparing  densities  in  one  dimension  as  a  function 
of  T.  It  is  certainly  much  taster  than  computing  in  two  dimensions  and  serves  as  an  effective 
rejection  filter.  Properly  used,  it  might  be  possible  to  do  all  or  nearly  all  calculations  with 
one-dimensional  densities.  One  can  extend  the  idea  and  define 

Q[t)  =  5(/)  -1-  S{t  +  T)  +  S{t  +  2T).  (14) 

The  density  associated  with  Q[i)  has  yet  a  different  structure.  One  conld  also  use  a  different 
7  ,  say  T\  instead  ot  'IT .  All  ot  thc.se  concepts  need  to  be  examined.  We  now  know  that  thev 
all  work  to  some  extent  aiul  the  goal  is  to  refine  and  optimize  their  use. 

In  order  to  illustrate  the  possibility  ot  discrimination  with  one-dimensional  densities,  we 
chose  two  signals  and  added  them  together  with  cc|ual  amplitudes.  In  Figure  11  the  solid 
line  is  the  density  of  the  summed  signal.  The  dashed  line  is  the  fit  obtained  as  u.sual  by 
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using  p{k)  =  pi[k)'p2{k).  The  curves  would  exactly  coincide  if  the  length  of  the  signal  were 
sufficiently  long.  The  difference  in  the  two  curves  provides  a  measure  of  what  is  the  ’‘best" 
that  one  can  do  if  one  tries  candidate  signals  to  make  up  the  received  signal.  In  this  case  the 
most  stringent  test  was  to  try  a  signal  with  the  same  power  spectrum  as  one  of  the  original 
signals,  but  with  randomizetl  Fourier  phases  so  that  it  is  noise.  The  result  is  shown  in  Figure 
12.  .Although  the  fit  looks  good,  they  are  to  be  compared  with  those  in  Figure  11.  There  is 
clearly  a  difference  for  all  values  of  V’  and  oik'  could  easily  infer  the  correct  signal  even  with 
this  one-dimensional  construction. 

The  one-dimensional  tlensities  are  just  simple  projections  of  two-dimensional  or  higher 
densities.  However  there  art.-  alternative  projections  that  are  also  simple.  For  instance,  for 
the  Fourier  transform  in  two  tlimensions,  p(k),  one  can  fix  the  magnitude  of  k  and  vary  the 
polar  angle.  We  have  done  some  limited  testing  of  this  idea  and  the  ability  to  classily  seemed 
comparable  to  the  one-dimensional  densities  described  immediately  above.  The  reason  for 
considering  these  alternatives  is  that  it  is  desirable  to  use  a  projection  that  has  considerable 
structure,  and  the  region  around  =  0  is  geiK'rally  pretty  structureless.  Ordinary  [projections 
have  to  use  that  region  whereas  this  projection  is  specifically  devised  to  avoid  it.  It  also  has 
the  feature  that  it  is  simple  to  calcula.t(‘  and  the  interpolations  required  because  of  the 
unknown  signal  strength  are  also  simphx 

It  is  possible  to  extend  the  definition  comparable  to  R{t)  or  Q{t)  to  two  dimensions. 
We  have  done  that  in  a  limited  way  and  the  comparisons  seem  to  work  better;  however,  not 
much  testing  has  yet  been  done. 

.An  important  new  tool  that  we  have  been  developing  is  the  use  of  an  array  of  receivers. 
We  use  the  array  first  to  get  the  direction  of  a  particular  source  and  then  to  analyze  the 
source.  We  made  a  model  of  an  array  typical  of  a  possible  underwater  array.  The  particular 
test  we  ran  was  for  one  hundred  elements,  spaced  a  few  meters  apart.  The  details  of  designing 
an  ap|)ropriate  array  are  pre.sented  in  the  next  section. 

Both  equal  spacing  and  logarithmic  spacing  were  tried.  Two  sources  a  long  distance 
away  with  an  angle  of  -ly  separating  them  were  simulated.  In  one  example  we  made  source 
A  ten  times  as  strong  as  B  (20db  separation).  We  were  easily  able  to  identify  B  with  a  single 
two-dimensional  comparison  or  a  few  one-dimensional  comparisons.  We  could  not  do  this 
without  the  array.  In  order  to  demonstrate  the  above  comments  we  show  the  results  with 
and  without  an  array.  Without  an  array  it  is  just  the  addition  of  two  signals  with  one  signal 
being  ten  times  larger  in  amplitude.  In  Figure  13a  we  show  the  density  the  sum  of  the  two 
signals, 

S{()  =  Si  +  lOS-i.  (15) 

In  Figure  Fib  we  show  the  density  obtained  by  the  usual  multiplication  of  pi  and  p2.  In 
Figure  13c  we  show  an  alternative  candidate  for  Si.  Clearly  it  is  difficult  for  the  eye  to  see 
the  difference  between  Figure  13b  and  Figure  13c.  It  is  also  numerically  diffiicult  if  only  one 
time  lag  is  used.  In  Figure  Fid  we  show  tlie  density  for  an  alternative  signal  for  S2-  Ihis 
time  we  can  clearly  see  the  dilference.  In  Figure  14  we  show  the  results  of  using  an  array. 
This  time  we  are  ch'aiiy  abh*  to  ideiitily  correctly  the  weak  signal.  Note  that  since  the  signal 
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is  broadband,  we  can’t  use  the  array  to  select  a  particular  frequency  in  a  given  direction  as 
one  might  do  with  harmonic  signals.  This  is  a  very  important  result  since  we  don’t  have  to 
search  for  pairs  of  signals,  etc.  but  only  need  to  match  a  single  signal,  as  the  other  received 
signals  were  greatly  reduced  in  strength. 

4.  OVERVIE’W  OF  BROAD  BAND  ACOUSTIC  ARRAY  DESIGN 

We  studied  the  problem  of  localizing  and  tracking  the  source  of  a  received  underwater 
acoustic  signal,  given  that  the  latter  is  more  comple.x  than  a  simple  CW  tone,  and  that  it 
has  propagated  from  source  to  receiver  through  a  non-linear  medium.  The  approach  we  had 
to  take  is  a  dual  one:  it  depends  on  how  well  one  can  represent  the  detected  signal.  Indeed, 
if  the  signal  produced  by  the  source  is  truly  chaotic,  i.e.  if  it  is  has  chaotic  dynamics,  then, 
while  we  have  an  approach  that  can  identify  the  dynamical  system  in  question,  one  cannot 
hope  to  retrieve  the  actual  time-series  ol  the  signal  with  any  kind  of  accuracy.  On  the  other 
hand,  if  the  signal  produced  by  the  source  is  merely  comple.x,  i.e.  if  it  is  more  complicated 
than  a  (linearly-produced)  time-harmonic  tone  yet  has  Hamiltonian  dynamics  or  at  least 
non-chaotic  ones,  then  one  can  indeed  hope  to  identity  a  template  time-series  associated 
to  it,  perhaps  parametrized  by  some  of  the  state  variables,  then  estimate  these  parameters 
along  with  the  other  state  variables  associated  to  the  source. 

In  either  case,  the  signal  produced  by  the  source  gets  convolved  by  the  Green’s  function 
of  the  medium  before  being  measured  by  the  receivers.  Because  the  source  is  moving  in  time, 
and  because  the  dependence  of  the  Green's  function  on  range  is  very  much  non-linear,  the 
received  signal  is  quite  comple.x  in  both  cases.  Yet  one  would  expect  that  il  this  complexity 
is  well-understood  and  optimally  used,  one  could  infer  from  it  the  location  ol  the  source  as 
a  function  of  time.  In  turn,  it  is  natural  to  expect  that  this  localization  and  tracking  should 
help  improve  the  detection  and  estimation  of  the  various  source  parameters. 

In  the  non-chaotic  case,  one  can  try  to  model  the  template  associated  to  the  source 
signal  as  accurately  as  possible;  indeed,  it  will  depend  on  the  source  coordinates  (which  vary 
with  time),  on  its  (time-varying)  velocity,  and,  possibly,  on  some  other  (non-time-varying) 
parameters  that  can  distinguish  one  underwater  source  from  another.  Once  the  dynamics 
of  these  state  variables  are  modeled  by  mathematical  e<[uations  (with  deterministic  and 
stochastic  components),  and  once  the  dependence  of  the  received  signal  on  these  variables  is 
also  modeled  mathematically  (i.e.,  essentially,  once  the  medium  is  satisfactorily  represented 
and  its  Green’s  function  computed),  we  reduce  the  problem  to  a  classical  optimal-filtering 
problem,  albeit  a  non-linear  one.  The  approach  we  take  in  this  case  is  to  derive  and  solve 
the  equations  for  the  conditional  density  function.  This  is  the  probability  density  function 
of  the  state  variables  describing  the  problem,  conditioned  on  the  signal  received.  Computing 
this  p.d.f.  is  essentially  the  equivalent  of  implementing  an  infinite-dimensional  Kalman  filter 
to  provide  the  best  estimate  of  the  state  variables  in  this  non-linear  problem. 

In  the  truly  chaotic  case,  one  cannot  expect  to  single  out  an  actual  signal  template.  \et 
the  approach  we  took  to  the  identification  problem  does  associate  a  characteristic  density 
function  to  our  source,  one  shared  by  all  the  diverging  actual  signals  that  could  be  produced 
by  this  single  dynamical  system  that  hevs  been  identified.  .As  this  density  function  scales 
linearly  with  amplitude,  we  can  also  estimate  the  instantaneous  (real)  amplitude  of  our 
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signal.  In  order  to  develop  a  tracking  method  in  this  case,  the  challenge  is  two-fold:  use 
plane-waves  (rather  than  the  detailed  Green’s  function)  to  make  the  best  use  of  what  little 
phase  information  one  has,  then  try  to  refine  this  angle-tracking  by  making  optimal  use  of 
the  sequence  of  estimates  of  the  instantaneous  real  amplitude  information  we  obtain  using 
the  density  function  that  characterized  the  signal.  Our  tracking  approach  in  the  chaotic 
situation  is  thus  somewhat  weaker  than  in  the  previous  case:  instead  ol  making  optimal  use 
of  the  complete  complex  Green's  function,  we  first  account  for  all  phase  information  using  a 
plain-wave  approximation,  then  account  for  the  varying  amplitude  using  the  real  propagation 
loss  function  (the  magnitude  of  the  Green’s  function). 

In  the  next  section,  we  describe  our  filtering  approach  in  the  chaotic  case  in  more 
detail,  and  present  simulated  results.  We  then  describe  our  approach  in  the  non-chaotic-yet- 
nonlinear  case  and  present  simulated  results  in  that  case. 

5.  BROAD  BAND  ARRAYS 

■As  explained  in  the  previous  section,  we  first  seek  to  extract  all  phase  information  from 
our  signal  before  processing  it  through  our  dynamical  .system  identifier.  lollowed  by  our 
instantaneous  amplitude  estimator.  I'lic  resulting  amplitude  estimate  is  then  used  along 
with  the  phase  information  as  tlu'  latest  data  for  our  tracking  algorithm.  We  now  describe 
the  individual  pieces  of  this  method. 

In  the  plane-waves  approximation,  the  best  we  can  hope  to  do  is  to  translate  the  phase 
information  into  a  bearing  angle.  The  problem  is  thus  to  design  an  ari'ay  whose  beam 
pattern  has  little  or  no  frequency  dependence.  In  particular,  we  must  make  sure  that  the 
main  lobe  width,  side  lobe  level,  and  distance  between  the  main  lobe  and  the  side  lobe  plateau 
(which  takes  the  place,  in  our  approach,  of  the  grating  sidelobes  which  appear  generically 
in  the  case  of  single  frequency  arrays),  as  well  as  the  peak  response,  all  remain  constant 
across  the  design  frequency  band.  To  do  this,  we  use  the  Poisson  summation  formula  and 
the  method  of  stationary  phase  to  relate  beam-pattern  properties  to  array  characteristics. 
Using  this  relation,  we  then  translate  all  the  listed  beam-pattern  requirements  into  practical 
requirements  on  hydrophone  location  and  array  amplitude  shading.  This  approach  is  similar 
to  the  one  proposed  by  Ishimaru  in  [3]  [4],  and  differs  in  that,  while  Ishimaru  et  al  only  verified 
that  a  particular  shading  function  gives  a  good  sidclobe  behavior,  we  use  the  approach 
to  systematically  derive  the  array  requirements  in  order  to  achieve  our  global  broadband 
objectives.  Our  generalization  of  his  method  to  the  systematic  broadband  case  here  follows 
the  same  line  as  our  generalization  of  his  method  to  time-domain  pulse  trains  ([2]). 

Let  us  begin  with  some  design  considerations.  The  problem  with  using  a  classical  uni¬ 
form  element  spacing  comes  from  the  two  facts  that 

•  in  order  to  reduce  the  side  lobe  levels  across  the  band,  the  spacing  would  have  to  be 
proportional  to  the  wavelength  corresponding  to  the  highest  design  frequency,  while 

•  in  order  to  keep  the  width  of  the  main  lobe  constant  over  the  design  band,  the  length 
of  the  array  would  have  to  be  proportional  to  the  wavelength  corresponding  to  the 
lowest  design  frequency. 
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In  other  words,  with  unil’orm  spacing,  the  total  number  of  elements  would  be  proportional 
to  the  ratio  of  highest  to  lowest  design  frequencies.  With  our  design,  it  turns  out  that  the 
number  of  elements  will  grow  like  the  logarithm  of  this  ratio,  hence  should  prove  substantially 
more  efficient. 

The  main  difficulty  in  using  unequal  spacings  is  that  the  beam  pattern  function  for  such 
an  array  is  unlikely  to  have  a  closerl  form  e.xpression.  This  makes  it  difficult  to  relate  array 
beam  pattern  properties  to  array  specihcatioiis  (i.e.  spacing  and  amplitude  shading).  That 
is  where  the  Poisson  summation  formula  comes  in. 

Suppose  we  have  iV  +  1  elements,  to  l)c  placed  along  tlie  .r-a.xis  at  positions  Xn  (0  < 
71  <  N),  and  to  be  shaded  by  the  (frecpiency-dependent)  weights  n’„(A),  A  a  wavelength 
within  the  band  of  interest.  As  a  function  /?  of  o  =  sin  (angle  of  the  plane  wave  arrival) 
—  sin  (array  beam  angle),  the  r(?sponse  of  the  array  will  be 

,v 

/?(a,A)  =  ^u.,.(A)e--^‘""’‘/'  (16) 

H=0 

There  are  many  ways  of  extending  .r„  and  (e„(A)  to  functions  of  a  real  argument  n  such  that 


•  both  functions  are  continuous  at  all  the  integers. 

•  X  is  twice  differentiable  and  strictly  increasing  on  the  interval  [0.  A  ]. 

•  w  vanishes  for  n  not  in  ]  — t.  .V  +  <[  (where  c  is  a  small  ])ositive  number  fixed  once  and 
for  all). 


Choose  such  extensions  x{n).  w(7uX).  We  can  then  use  the  Poisson  summation  formula  on 
R  to  get 


y?(o,A) 


,c(7i,  dll 


(17) 


For  simplicity,  let  us  assume  that  the  element  spacings  are  increasing  (with  n).  Mathemat¬ 
ically.  we  are  assuming  that  the  derivative  x'  is  an  increasing  function,  i.e.  that  the  second 
derivative  x''{n)  >  0  for  n  in  the  interval  [0,  iV].  Following  [3]  [4],  let  us  use  the  method  of 


stationary  phase  to  study  the  expression  (17).  This  method  predicts  that,  asymptotically, 
(i.e.  for  A  relatively  small),  the  behavior  of  R  for  a  close  to  0  will  be  adequately  described 
by  that  term  in  the  right  hand  side  of  (17)  which  corresponds  to  rn  =  0  (i.e.  the  mainlobe 
is  described  by  the  m  =  0-term),  wliile  the  grating  sidelobe  will  be  described  by  the 
behavior  of  the  integrals  near  that  n  =  n„,  which  gives  a  stationary  phase  in  the  integrand, 
i.e.  for  which 


(IS) 


Since  x”  was  assumed  |)ositiv('.  x'  is  strictiy  increasing,  so  (IS)  will  have  at  most  one  solution 
for  each  in. 
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Let  us  begin  by  examining  the  mainlobe-to-sidclobe  separation  problem.  What  is  the 
smallest  positive  value  ot  a  for  which  (IS)  has  a  solution,  i.e.  where  does  the  first  sidelobe 
appear?  For  a  sufficiently  near  0,  the  right  hand  side  of  (IS)  is  very  large,  outside  the  finite 
range  ,r'([0,  .Vj)  of  .r'.  A  s  a  increases,  the  right  hand  side  will  decrease  accordingly.  The  first 
solution  thus  occurs  when  a  =  /7i.A/;r'( .V).  Hence,  at  the  lowest  design  frequency,  where  we 
can  and  will  make  the  reasonable  assumption  that  the  whole  array  is  used,  the  first  (?7i  =  1) 
sidelobe  appears  when  a  =  oq  =  Au/.r'(//),  where  Ay  is  the  wavelength  corresponding  to  the 
lower  design  freciuency,  and  .r^(A  )  is,  by  hypothesis,  the  largest  element  spacing.  In  order 
to  meet  our  requirements,  we  would  like  to  keep  this  main-to-side-lobe  separation  constant 
across  the  design  band.  The  obstacle  is  that,  if  A  <  Aq,  (IS)  will  have  solutions  with  a  <  qq. 
The  only  wa\'  to  remo\’e  this  hurdle  is  to  make  sure  that  the  corresponding  shading  weight 
equals  zero,  for  all  these  potential  .solutions.  So  we  have  to  set 

fc(n,A)  =  0  for  all /I  >  x'  '(x-'(.V)A/Ao)  (19) 

and  keep  it  positive  otherwise.  This  equation  has  a  nice  interpretation:  it  says  that  at 
any  given  frequency  in  the  design  band,  elements  spaced  farther  than  .r'(iV)/Ao  wavelengths 
apart  should  not  be  used. 

Let  us  now  examine  the  sidelobe  level  problem.  .As  in  [l],  (17)  implies  that  the  level  of 
the  m‘^  grating  sidelobe  is  approximately 

I  /  a’(n,A)e--^'‘^“'*“>/'-"‘'‘>dn|*  (20) 

J  —  >0 

which,  by  stationary  phase,  is  it.self  approximateiy  equal  to  Aaq??,,,,  A)’/(Qi'"(7?,n)),  with 
n-m  —  x'  *(77iA/a)  as  in  (IS).  That  is,  the  level  of  the  sidelobe  is  approximately 

1  „  ,fC-(77,„,A)2 

— -i-(n.m) — - —  21 

rn  .r"(/u) 

.As  we  determined  in  the  previous  paragraph,  the  first  (and  largest)  sidelobe  occurs  at  q  = 
A/.r'(A'’)  >  c7o  =  Ao/x'(.\').  i.e.  for  a  in  the  interval  [Ao/x'(A'’)  ,  2],  for  which  ni  must 
correspondingly  lie  in  the  interval  [ir'  *(A/2),.r'  *(:c'(7i)A/Ao)].  For  this  sidelobe  level  to  be 
independent  of  A,  we  have  to  require  that 

iv{n,\)^  =  constant— <  n  <  x'~^{x'(N)X/2)  (22) 

.r'(77.) 

over  all  A  in  the  design  band. 

Let  us  now  turn  to  the  mainlobe  level  problem.  .As  in  [1][3],  the  array  response  for  a 
near  0  is  approximately  given  by  the  0*^'  term  in  (17),  namely 

/  U7(77,A)c""'""^'‘^/'f/77  (23) 

J  — OO 

Choosing  t  =  x{n)  as  a  new  variable  as  in  [2],  this  can  be  rewritten  as 


»’(.!■  '(/),A)e 


X'(2-1(0) 
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where  we  have  made  the  reasonable  assumption  that  the  array  extends  from  0  =  x(0)  to 
x{N)  =  lo-  There  are  several  ways  to  make  sure  that  this  integral  is  independent  of  A,  at 
least  for  a  small.  One  simple  way  is  to  write  L{X)  for  the  length  of  that  part  of  the  array 
used  at  wavelengtii  A.  i.e.  to  assume  that  Lv(n,\)  vanishes  exactly  when  n  >  a.'~'(Z,(A)),  so 
that  the  integral  can  be  rewritten  as 


rLiX) 

Jo 


iu(x 


.-1 


x'(x-'{t)) 


a/.\ 


(25) 


then  to  assume  that  the  integrand  is  as  spatially  uniform  as  possible,  i.e.  depends  only  on 
A.  VV'e  thus  assume  that  »’(.?•“*(/),  A)  /  ;r'(.r~‘ (/))  =  /(A),  and  there  remains  to  determine 
the  functions  /  and  L.  VVe  do  this  by  v<Mifyiug  when,  with  llu^sc  functions,  the  integral  is 
indeed  independent  of  A: 


/(A)e-’-‘‘'/V/  =  /(A) 


g2,T,aL(,\)/A  _  I 

'Ixiaj  \ 


(26) 


which  is  indeed  independent  of  A  only  if  L(A)  =  LoX/Xq  (we  have  already  assumed  that 
the  entire  length  of  the  array  Lq  is  used  at  tlio  lowest  design  wavelength  Ao),  and  if  /(A)  = 
constant' j X.  We  summarize  these  conclusions: 


w[n,X) 


constant' 


An) 

X 


if  72  <  x-“^(ZoA/Ao), 


and  should  vanish  otherwise. 


(27) 


What  should  the  overall  design  be  in  light  of  all  this?  The  first  halves  of  conditions  (19) 
and  (27)  specify  the  support  of  the  shading  function  w.  To  reconcile  them,  we  must  succeed 
in  construction  x{n)  such  that 


Ao  Ao 


(2S) 


Using  the  new  variable  A'  =  LqX/Xq,  this  ecpiation  becomes  .i-“'(A')  =  x'~^{x'{N)X' Lq). 
If  we  now  write  u  for  .x-“'(A')  (so  that  A'  becomes  x{i/))  then  apply  x'  to  both  sides  ot  the 
equation,  we  obtain 


dx  _  x'(N)  _ 
du  Lq 


(29) 


a  differential  equation  which  is  solved  by  x{u)  =  ^  which  .4  is  a  constant  that 

can  be  determined  by  the  condition  x{N)  =  Lq.  Putting  it  all  together.  The  solution  is 


.r„  =  Lqc  ‘■0 


(30) 


This  equation  must  hold  at  least  for  all  n  (=  x~^(ZoA/Ao))  that  are  within  the  interval 
[Ai,  Ao]  where  Ai  is  the  wavelength  corresponding  to  the  highest  design  frequency.  Also  note 
for  future  use  that  (30)  implies  that,  over  its  range  of  validity.  .r(.T'~*((/))  =  LQd/x'{N). 


With  formula  (30)  for  x^,  condition  (22)  governing  the  sidelobe  levels  becomes 

((;(;?,  A)'  =  constant  I / Lo  (31) 

for  all  71  in  the  range  2''“^(/\/2)  <n  <  x'~\x'{N)\/ Xo)  and  all  A  in  [Ai,  Aq].  It  is  important 
to  note  that  in  order  to  validate  the  use  we  just  made  of  formula  (30),  we  must  require 
that  that  formula  hold  for  all  n  in  the  interval  [x'  '(A/2),x'  '(x'(.V)A/Ao)]  over  all  A  in  the 
design  band.  Since  x'  is  increasing,  this  means  that  (30)  must  hold  for  all  7Z  in  the  interval 
[x'  ^ (Ai/2),  iV].  As  to  our  amplitude  shading  function  w,  we  ma\’  as  well  normalize  it  by 
setting  the  constant  in  (31)  equal  to  Lo,  so  that  w{n,X)  —  1  over  the  interval  above. 

To  specify  over  the  entire  domain  of  interest,  we  turn  to  the  last  condition  unaccounted 
for  yet,  namely  the  second  half  of  (27),  namely  that  w{n,X)  =  constant'  x'{n)/ \  if  < 
x“^(ZoA/Ao).  Unfortunately,  to  maintain  a  low  sidelobe  level,  we  have  already  adopted  the 
condition  w{7i,X)  =  1  for  .r~^(A/2)  <  n.  <  x“‘(/.oA/Ao).  So  we  will  simply  choose  the 
constant  in  (27)  in  such  as  way  as  to  make  te  continuous  at  the  transition  point  n  =  x“'(A/2), 
i.e.  we  set  constant'  =  2,  and  set 

7o{n.X)  =  ^.x'(?r)  for -n  <  x'  '(  — )  (32) 


To  obtain  a  complete  design,  we  still  need  to  specify  ,i'„  for  all  =  1  ....  A'.  So  far  we 
know  that  (30)  must  hold  for  all  n  in  tiie  intersection  I  of  the  intervals  [x“’(Io'^i/'-\o))  A] 
and  [x'~^ (Ai/2),  A^].  Therefore  we  still  need  to  determine  .'r'(.'V)  and  the  values  of  x,i  for  n 
outside  I. 


.Assuming,  without  loss  of  generality,  that  x'(A^)  >  Ao/2  (i.e.  that  the  largest  hy¬ 
drophone  spacing  is  at  least  as  big  as  half  the  largest  wavelength  in  the  problem),  our  range 
I  is  in  fact  the  interval  I  =  [x'~*(Ai/2),  A^].  Thus,  (30)  specifies  the  element  positions  from 
the  widest-spaced  down  to  the  ones  spaced  at  half  the  shortest  design  wavelength.  The  re¬ 
maining  elements  (whose  positions  are  yet  to  be  specified)  must  therefore  be  spaced  at  least 
as  closely  as  half  the  wavelength  at  any  frequency  in  the  design  band.  For  efficiency,  it  is 
only  natural  to  minimize  the  number  of  additional  elements  available  to  fill  out  the  array, 
and  therefore  require  that  they  be  spaced  no  closer  than  the  half-minimum-wavelength.  So 
we  set 


-•^1  r  ^  -1/  Lo  Ai  Lo  ,  x'{N) 


(33) 


In  order  that  x  remain  continuous,  we  must  finally  reconcile  this  last  equation  with  equation 
(30)  at  the  intersection  point  of  their  ranges.  This  reduces  to  imposing  the  final  additional 
condition 


-^1 

~ 


Lq 

2  x'{N) 


(34) 


Putting  all  thc.se  requirements  together,  we  end  up  with  the  following  equations  for 
the  element  positions  and  the  shading  coelTicients  (here  we  call  d  the  largest  spacing  in  the 
array): 

d  >  Ao/2  (35) 


IS 


1  +  log(2c//A,  ^ 

.i-„  =  /;Ai/2ror;).  <  N  -  ^  log(2J/Ai) 

a 

Xn  =  lor  n>N-^  log(2f//Ai) 

a 

Wn(\)  —  2x'{n)/X  for  ?2  <  iV - y  log{2d/\) 

a 

WniX)  =  1  for  iV  -  ^  log(22//A)  <  n  <  log(Ao/A) 

d  d 

"’n(A)  =  0  otherwise. 

Practically,  this  implies  that  the  filtering  we  apply  to  the  Fourier  coefficients  I'nif)  at 
freciueucy  /  of  the  signal  at  the  /jf''  hy(lroj)hone  is 

(42) 
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with  as  above,  /l„(/)  -  \i/C  if  n  <  N  -  ^  log(2(//Ai ),  /!„(/)  =  (2c//C)e("-''')‘'^» 
if  N  —  •^log(2J/Ai)  <  N  —  log(2(///C') ,  and  0  otherwise,  while  Bn{J)  =  1  if  ds  — 

df-  log(2f///C')  <  n  <  N  —  log(Ao//(") ,  and  0  otherwise,  and  ,3  is  the  sine  of  the  bearing 

angle  being  scanned. 

As  the  simulations  that  wc  have  implemented  show,  this  method  should  prove  very 
useful  in  concentrating  the  energy  of  a  signal  spread  over  a  wide  band  and  thus  enhancing 
its  detectability.  Indeed,  the  method  singles  out  the  optimal  set  of  time  delays  (defined  by 
the  angle  of  peak  gain)  to  use  in  adding  up  coherently  the  signals  arriving  at  the  different 
elements  of  our  array. 

VVe  have  previously  shown  the  benefits  of  an  array  in  signal  classification  (see  figures 
13  and  14).  In  the  next  example  we  show  how  one  can  determine  the  direction  of  a  broad 
band  signal,  even  when  a  stronger  signal  is  present.  In  Figure  15  we  show  an  example  with 
two  sources,  one  at  —45°  and  one  at  +30°.  The  array  does  an  excellent  job  of  locating  the 
direction.  Knowledge  of  the  direction  is  important  for  optimum  utilization  of  the  array  for 
classifying  signals.  Indeed,  the  method  singles  out  the  optimal  set  of  time  delays  (defined  by 
the  angle  of  peak  gain)  to  use  in  adding  up  coherently  the  signals  arriving  at  the  different 
elements  of  our  array. 

6.  NONLINEAR  TRACKING 

In  this  case,  we  are  assuming  that  a  “parametrized”  signal  template  time-series  has  been 
identified.  By  “parametrized,”  wc  mean  that  this  template,  which  we  have  as  a  model  for 
the  signal  that  we  are  trying  to  detect  and  localize  in  tlu2  received  time  series,  may  depend 
on  several  parameters  such  as  target  aspect,  speed,  range,  etc.  VVe  also  assume  that,  while 
the  dynamics  of  this  template  signal  are  not  necessarily  linear,  they  are  not  chaotic:  in  other 
words,  we  can  expect  to  estimate  in  a  robust  fasliion  the  values  of  the  parameters  on  which 


(36) 

(37) 

(38) 

(39) 

(40) 

(41) 
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it  depends,  and,  indeed,  locate  the  exact  version  of  onr  signal  within  the  received  time  series. 
For  convenience,  let  us  group  all  our  parameters  into  a  “state”  vector  X,  and  assume  that 
we  are  looking  for  the  signal  S{t\X)  within  the  received  time  series  r{t). 

So  far,  we  have  not  accounted  for  the  (highly  non-linear)  propagation  effects  through 
the  ocean  lens.  For  simplicity,  we  will  assuiiu'  for  now  that  the  parameter  vector  .V  consists 
only  of  the  coordinates  of  the  source  producing  the  time  series.  Our  problem  in  this  case  is 
to  look  for  each  of  the  n  signals 

Sn[l;.\)  =  j  C'(.\ ,  Iq]  X,i,  t)  S{to:  X)  dto  (4-3) 


within  the  corresponding  signal  received  at  the  hydrophone,  where  G  is  the  Green's 
function  for  the  medium,  and  .V„  represents  the  coordinates  of  the  array  element.  Rather 
than  try  to  detect  5  (and  estimate  its  parameters  X)  based  on  a  single  iook.”  an  approach 
whose  success  is  assured  only  when  one  assumes  an  unrealistically  high  signal-to-noise  ratio, 
we  propose  to  combine  several  looks  together  in  order  to  make  our  detection.  In  short,  we 
set  out  to  “track”  S  (and  A  )  in  order  to  d('terniine  if  it  is  indeed  there,  i.e.  for  what  value 
of  A"  ,  if  any,  can  we  identify  our  template  within  the  received  time  series. 

The  optimal  filtering  solution  to  this  problem  is  described  by  the  formalism  of  the  Zakai 
equation  ([6'j).  The  solution  to  this  cciuation  is  essentially  the  (infinite-dimensional)  non¬ 
linear  equivalent  of  the  linear  Kalman  filter.  To  set  up  the  Zakai  ecpiation.  and  then  try  to  find 
its  solution,  one  starts  by  describing  the  exi)cctcd  time  evolution  of  the  vector  one  is  trying 
to  estimate,  .A  in  our  case.  One  then  describees  the  dependence  of  the  data  received,  the  I'Gs 
in  our  case,  on  this  state  vector  -  that  is  achieved  by  our  equation  (43).  once  we  can  compute 
the  medium's  Green’s  function.  For  simplicity,  we  shall  assume  very  simple  dynamics  for 
our  variables,  and  proceed  to  write  down  directly  the  solution  to  the  Zakai  equation:  first, 
we  change  the  assumption  about  A'  a  little,  and  assume  that  our  state  vector  consists  of  the 
coordinates  of  the  source  as  well  as  its  velocity  vector,  i.e.  X  =  (.r*’*,  i’^^*), 

then  we  assume  that  the  dynamics  of  the  velocity  vector  have  a  deterministic  part  that 
assumes  that  the  velocity  is  piecewise  constant,  and  a  stochastic  part  that  assumes  that 
the  discrete  velocity  changes  occur  at  Poisson-distributed  times  with  rate  /i.  This  measn 
that  our  source  moves  in  straight  line  segments,  changing  its  course  at  a  sequence  of  times 
{ii,  <2,  <3...}  such  that  the  quantities  txJ.-i  —  fiAa  ~  h-,  ■■■  Poisson  distributed  with  mean 
l//i,  and  such  that  the  new  velocities  at  each  course  change  are  distributed  according  to  a 
pre-specified  density  function,  call  it  g.  The  formalism  of  the  Zakai  equation  allows  us  to 
determine  the  conditional  density  function  p(TV’ given  /’^(r)  for  all  r  <  t),  the  function  of 
{t.Y)  which  describes  the  likelihood  that  A’  at  time  t  be  close  to  the  value  T,  given  all  the 
past  received  observations  r„{T). 

Under  these  assumptions,  we  have  obtained  a  recursive  algorithm  for  updating  p  as  the 
data  is  received.  Indeed,  if  the  noise  in  the  received  signal  is  0-mean  additive  Gaussina  noise 
with  r.m.s.  level  cr,  p  can  be  computed  using 


p(U.c,u) 


p{to,x  -  (t  -  tojV.VjC  ''0 


(44) 
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+/'  /  [  /  pito,:t-[t-to)v,v  )g{v-v  )dr\c  '''0 


where  u„(i’)  ==  rn{T)dT  is  the  sum  uf  ;ill  ihita  received  at,  the  /r‘**  element  up  to  time  .s. 


We  set  out  to  study  the  moments  of  p  and  tlieir  evolution  in  time  in  order  to  determine 
the  kind  of  localization  one  can  expect  to  achieve  efficiently.  The  preliminary  simulations  we 
had  conducted  started  by  a.ssuminga  nK'dium  with  given  visco-elastic  propagation  properties, 
then  computing  the  associated  Green's  function  for  a  pure  tone,  thus  enabling  us  to  compute 
the  templates  S'n  foi'  '"my  array  conliguratiou. 


The  first  new  aspect  of  the  problem  we  (hxided  to  study  was  the  robustness  of  this  ap- 
proacli  when  one  does  not  necessarily  know  the  exact  parameters  governing  the  propagation 
in  the  medium  at  hand.  Indeed,  oiu'  cannot  expect  to  know  the  compressional  and  shear 
velocity  profiles  with  all  their  fluctuations  throughout  the  medium,  nor  can  one  expect  to 
compute  an  exact  Green’s  function,  even  iftlie  propagation  parameters  were  known  e.xactly. 
In  particular,  starting  with  a  uniform  /i  over  a  realistic  domain  for  A’,  the  first  question  is 
to  determine  how  soon  the  variance  of  p.  i.<\  the  mean  sriuare  nneertainty  in  the  estimate 
of  ,Y,  would  converge  to  an  acceptably  small  value.  We  started  by  doing  this  tor  a  pure 
tone  (the  case  of  a  linear  signal),  in  a  two-layer  medium  allowing  only  compressional  wave 
lu’opagation.  The  discouraging  laxult  is  that  the  nu-thod  fails  if  the  relative  discrei)ancy 
between  the  modeled  velocity  and  the  actual  propagation  speerl  is  greater  than  O.  f  u.  if  the 
relative  difference  between  tlie  (riod<“l<>d  and  actual  !;{}<'{'  drpll!  is  greater  tluin  .Jit.  or  if  tite 
relative  error  in  the  density  of  the  half-infinite  layer  is  greater  than  109c.  .Act  only  rloes  the 
second  moment  fail  to  converge,  but  the  algorithm  actually  eiuls  up  "zeroing"  p  because  the 
increasingly  negative  exponents  in  (-14)  impose  a  uniform  rapid  decay  rate  on  the  density 
function.  In  essence,  the  "ambiguity  function”  for  this  particular  approach  is  too  narrowly 
supported,  and  this  model  is  too  sensitive  to  randomness  in  the  medium. 


The  conclusion  is  that  for  this  approach  to  work,  one  must  incorporate  into  one's  vector 
of  “state  variables”  the  parameters  that  describe  the  medium,  and  one  must  incorporate  the 
effect  of  their  fluctuations  on  the  signal  as  it  pro|)agates  from  the  source  to  the  receivers. 
This  is  considerably  more  complicated  than  what  we  set  out  to  study:  such  an  endeavor 
would  lie  beyond  the  scope  of  this  effort. 
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Figure  1.  Time  traces  of  four  of  the  signals  used  in  this  demonstration.  Signal  3  is  random 
gaussian  noise  with  the  same  spectrum  as  one  of  the  other  signals. 
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Figure  2.  Time  traces  of  four  of  the  signals  used  in  this  demonstration.  Signal  7  is  random 
gaussian  noise  with  the  same  spectrum  as  one  of  the  other  signals. 
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Figure  3.  The  corresponding  power  spectra  for  the  signals  shown  in  Figure  1 


Figure  4.  The  corresponding  power  spectra  for  the  signals  shown  in  Figure  2 
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Figure  5.  The  density,  p{x,  y),  for  a  particular  signal  is  shown  in  Figure  5.  The  density  from 
another  sample  of  the  signal  is  shown  in  5c,  and  the  density  from  a  different  signal 
is  shown  in  5b.  These  latter  two  densities  are  examples  of  proposed  candidates  taken 
from  a  library  of  possible  densities. 
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Figure  6.  (a)  is  the  power  spectrum  of  the  signal  shown  in  6b;  (b)  is  the  sum  of  signals  4  +  6 
with  the  relative  strength  parameters  a  and  b  given  in  the  text;  (c)  Different  pair  of  signals 
that  was  being  compared  with  signals  4  +  6;  (d)  is  a  figure  of  merit  for  how  good  the  fit 
is,  smaller  x^  being  better.  The  horizontal  axis  is  a  relative  amount  of  the  two  signals  being 
combined.  The  two  curves  that  nearly  touch  are  two  samples  of  signals  4  +  6.  The  upper 
curve  is  for  signal  6  +  another  unshown  signal,  and  the  x^  for  signals  1  +  5  would  be  off  the 
graph. 
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Figure  7.  Figure  7a  is  the  density,  p{x,y),  for  S  =  .375'i  +  .635'2  with  Si  and  S2  being  two  of  the 
signals  described  above.  Figure  7b  is  the  density  constructed  from  p{k)  =  pi(.37/c)p2(-63k) 
and  would  be  identical  to  Figure  7a  for  an  infinitely  long  time  series.  Figures  7c  and  7d  are 
constructed  similarly  to  Figure  7b,  but  the  signal  used  for  Si  is  incorrect. 
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Figure  8.  One-dimensional  density  for  one  of  the  signals  with  Ssum  (0  =  S{t)  +  S{t  +  T)  and 
various  time  lags,  T.  The  vertical  axis  is  the  density  and  the  horizontal  axis  is  the  signal 
amplitude. 


Figure  9.  Same  as  Figure  8  for  another  signal. 
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Figure  10.  The  signal  has  the  same  spectrum  as  the  signal  in  Figure  9,  but  the  Fourier  phases  were 
randomized  to  produce  a  noise  signal.  The  single  figure  on  the  last  row  is  for  S{t)  =  cos(t). 
The  densities  for  cos(t)  are  independent  of  the  time  lag,  T.  The  other  densities  have  very 
little  dependence  on  T. 
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Figure  11.  The  solid  line  is  the  density  constructed  from  S{t)  =  Si{t)  +  S2{t)  where  and 
S2  are  two  of  our  standard  signals.  The  dashed  line  is  the  fit  obtained  by  calculating  the 
convolution  of  the  densities  for  and  S2-  The  difference  is  due  to  the  limited  length  of  the 
signal. 


Figure  12.  The  solid  line  is  the  same  as  that  in  Figure  11.  The  dashed  line  comes  from  replacing 
Si  by  another  signal  (see  text)  and  then  calculating  the  convolution  of  its  density  with  /?2> 
the  density  for  signal  52- 
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Figure  13.  The  signal  is  S{t)  =  Si{t)  +  iOd'2(<)  and  Figure  13a  shows  the  density  for  S.  Figure 
13b  is  the  density  obtained  from  p  «=  Pi  To  obtain  Figure  13c,  we  construct 

a  trial  p  from  a  signal  other  than  5,  and  then  Fourier  transform  it.  Figure  13d  is  similar 
e.xcept  that  a  signal  other  than  Si  is  used.  We  see  that  with  a  signal  this  weak  w'e  can  not 
identify  it  with  a  simple  receiver  and  a  single  time  lag,  although  the  strong  signal  is  readily 
identifiable. 
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Figure  14.  Two  signals  were  directed  at  an  array  with  an  angle  of  45°  between  them.  Signal  1 
had  an  amplitude  10  times  signal  2,  and  the  array  had  100  elements.  The  array  was  pointed 
at  the  weak  signal.  The  two-dimensional  density  of  the  signal  is  shown  in  Figure  14b.  Two 
sample  hypotheses  are  shown  in  Figures  2  and  3.  These  two  signals  happen  to  have  identical 
power  spectra.  Clearly  Figure  14c  can  be  easily  selected.  Without  the  array  and  with  only 
one  time  lag  used,  we  could  not  have  identified  the  correct  component. 


